A = [0.75,0.5,0.25;0.125,0.25,0.5;0.125,0.25,0.25];
p0 = [0.5,0.25,0.25]';
tp = [];
for i = 1:100
    p0 = A*p0;
    tp(i,1:3) = p0;
end

[P,D] = eig(A)1
